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In this work we study ultracold Fermions confined in a two-dimensional optical lattice and we explore the 
Mott-insulator transition with the Fermi-Hubbard model. On the basis of a mean-field approach, we study the 
phase diagrams in the presence of a harmonic trapping potential. Local Mott-insulator phases are shown to be 
generally situated in the center of the trap and correspond to a vanishing variance of the local density. We then 
study the effects induced by rotation on the Mott-insulator phase transition. In particular, we show that the phase 
boundary reproduces the edge of the Hofstadter butterfly. 
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I. INTRODUCTION 

Today great efforts are devoted to the investigation of 
fundamental phenomena using ultracold atoms trapped in 
optical lattices d |2l. One of the most outstanding break- 
throughs in the field concerns the experimental realization of 
the bosonic Mott-insulator transition |3 |. The transition from 
the superfluid to the Mott-insulator is performed in ultracold 
bosonic optical lattices by increasing the ratio of the on-site 
interaction to the tunneling amplitude IK O. When the 
Mott-insulator phase is reached, the system is characterized 
by a commensurate filling of the lattice. As predicted by 
Jaksch et ah , this quantum phase transition can be achieved 
by varying the depth of the optical potential (61 . 

The fermionic Mott metal-insulator transition (MMIT) has 
been originally predicted for electronic systems (71 . Recently 
it has been suggested that a fermionic Mott-insulator phase 
should be observed in one-dimensional optical lattices |8, 9l. 
Rigol et ah have indeed shown, using quantum Monte Carlo 
simulations, that confined Fermions interacting through a 
repulsive interaction may undergo an incompressible state 
characterized by full occupancy of the sites. When the 
confining harmonic potential is taken into account, the 
Mott-insulator phase is generally restricted to a finite domain 
at the center of the trap surrounded by compressible metallic 
phases | 8 1. Liu et al. have put forward the fact that collective 
oscillations of the atomic mass density constitute efficient 
signatures of the MMIT f9l. This result, based on the 
Bethe ansatz solution of the one-dimensional Fermi-Hubbard 
model, is currently motivating experimentalists to investigate 
the fermionic Mott-insulator phase in optical lattices 1 10','TT|. 
Very recently, the MMIT has been theoretically investigated 
for a fermionic system with interactions near the unitary limit 
(T2l [T3l . The competition between the BCS-superfluid and 
the insulator phases in strongly interacting Fermions systems 
has also been studied in Ref.| 14|. Finally the MMIT has been 
investigated in the context of three-dimensional optical traps 

ca. 

In this context, it is interesting to study the behavior of 
the MMIT when the system rotates (161 El. or equally, 
when a light-induced effective "magnetic" field is produced 



in the optical lattice (El [HI [20l ED [22 EJI- Such artificial 
fields, which are currently realized in laboratories in order 
to simulate the dynamics of charged particles in a magnetic 
field, open the door for profound explorations of condensed 
matter physics (23l . In particular, optical lattices subjected to 
effective "magnetic" fields constitute versatile systems which 
allows the investigation of vortex-physics L24l llgl [26[ [27l l28l 
and the quantum Hall effects (2Q1 [29l [30l [SH [32. In 
the context of Rose gases, the superfluid-insulator phase 
transition is known to be affected by effective "magnetic" 
fields (28l [33] |34l. Quite surprisingly, the phase boundary 
for this strongly-correlated system highly depends on the 
single-particle energy spectrum. The latter, originally studied 
in the context of electrons subjected to high magnetic fields 
by Hofstadter and Wannier, is fractal and its representation as 
a function of the effective "magnetic flux" adopts the shape 
of an intriguing butterfly, the so-called Hofstadter butterfly 
(35l [3611 . It has been proved recently that the phase boundary 
for the superfluid-insulator transition is directly related to 
the Hofstadter butterfly's edge [W, "34]. In a very different 
context, this highly irregular curve is also known to represent 
the phase boundary for the normal- superconducting phase 
transition in superconducting networks (37[r38ii39i . 



In this work, we study the MMIT for a Fermi gas con- 
fined in a two-dimensional optical lattice. Contrary to the one- 
dimensional case, one has to treat the Fermi-Hubbard model 
on the basis of a mean-field approach. We first study the phase 
diagram in terms of the hopping amplitude and the chemi- 
cal potential. We show that Mott-insulator phases, character- 
ized by integer fillings and vanishing variance of the density, 
should be observed in two-dimensional fermionic optical lat- 
tices. In the presence of a harmonic trapping potential, we 
show that local Mott-insulator phases are generally present 
at the center of the trap, while the average variance is non- 
zero. We then study the system in the presence of an effective 
"magnetic" field. In particular, we study the phase boundary 
in terms of the effective "magnetic" flux and compare it to the 
single-particle spectrum. 
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II. MEAN FIELD THEORY 

We consider the case of a trapped fermionic gas in the pres- 
ence of an effective "magnetic" field. We assume that the opti- 
cal potential created by the lasers is sufficiently strong in order 
to apply a tight-binding approximation. The effective "mag- 
netic" field B = Biz, which is supposed to be created in the 
system with lasers ifTSI or by rotation |17ll25l, is character- 
ized by the parameter ^. The latter is the number of effective 
"magnetic" flux quanta per unit cell. For a rotating system, 
this parameter is given by <l> = 2Mo?Q./h, where Vt is the 
angular velocity, M is the particle's mass, a is the lattice con- 
stant and h = 27rh is Planck's constant (TT). We can treat the 
many -body problem by considering the Landau gauge ifTSl, 
A = (0,5x,0), in which the Fermi-Hubbard Hamiltonian 
reads 



^ frn-\-l,n,afm,n,a 



m,n 'm,n,a 



(1) 



where fm,n,a (fm n cr) ^hc fermionic annihilation (creation) 
operator of a particle on the lattice site {m^n) with spin 
^{=1^)^ and where the local density operator per spin is 
denoted by fim.n.a- The fermionic operators satisfy the an- 
ticommutation rule = ^{i — i')- The parame- 

ter t is the nearest-neighbor tunneling amplitude, U is the 
repulsive interaction (U > 0), /i is the chemical potential 
and Vharm is a trapping potential |2|. In recent experiments 
involving trapped ^^K, the optical lattice has lattice depth 
Vq = 5Eji, where = h? /2M\^ is the recoil energy, 
a = A/2 = 532 nm is the lattice constant, and the tunnel- 
ing amplitude t ~ x 290 Hz ifTTIl . The tight-binding regime 
is reached for Vq/Er > 2.23 flUl . which is the case for these 
experiments. In the following, we set the lattice constant a 
and the particle's mass M to unity and work in units where 
h=l. 

In the Landau gauge, the magnetic phase 0{m) is given by 



(m,n+l) 



Ay{m)dy = 27r^m, 



(2) 



where $ = 5/27r is the "magnetic" flux quanta per unit cell 
and where the integral is performed along the link connect- 
ing the lattice sites (m,n) and (m, n + 1) [ IS'l. We point 
out that a similar Fermi-Hubbard Hamiltonian can be obtained 
with the symmetric gauge A = x x ft and can be found in 
Refs. |[T7l[25l . In this work, the effective "magnetic" flux $ 
is supposed to vary between [0, 1]. This "high magnetic field 
regime" can easily be reached with light-induced gauge fields 
ifTSl . In the context of rotating optical lattices, this regime is 
reached for 1] 1 kHz. 

In the following we consider that the "magnetic" flux is ratio- 
nal, namely ^ = p/q where p and q are integers. In this case. 




Figure 1: (Color online) Hofstadter butterfly $ = ^(E): single- 
particle spectrum for non-interacting Fermions in a rotating optical 
square lattice. The butterfly's edge is colored in red. 



the single-particle Schrodinger equation (for U = Vharm = 0) 
yields the well-known Harper equation 1351 . This difference 
equation is expressed as 



ikx 



m+1 " 



ikx 



- 2 cos(27r$m + ky) ipm = Etjjj^ 



(3) 

where ipm is a g'-periodic wave function, k is the wave 
vector and E is the single-particle energy. The wave vector 
belongs to the magnetic Brillouin zone, a two-torus defined 
SiS kx G [0,^] and ky G [0,27r]. The energy spectrum 
associated to Eq. ([3j has a band structure, consisting of q 
subbands, which has been extensively studied in the literature 
1351 [36l l40l . The representation of the energy as a function 
of the "magnetic" flux <l> leads to a well known fractal, the 
so-called Hofstadter butterfly (Fig[T]). 



In order to study the interacting system Eq. ([T]), we approxi- 
mate the hopping term on the basis of the decoupling formula 
for the product of two field operators Il5l l24ll33l , 



fm,n,afm-\-l,n,a — {fm,n,a)fm-\-l,n,a ~^ fm,n,a{fm-\-l,n,a) 

{fm,n,a){fm-\-l,n,a)' (4) 



This approximation holds for small fluctuations defined as 

^m,n,a = {f^,n,a) " {fm,n,a)'^- The mean- field Hamiltonian 
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HuF then takes the form 



^ ^ ^ ( /m+l,n,cr^m,n,cr ^ ^m+l,n,cr/m,n,cr 
m,n,a ^ 

e ^ ^ (^m,n+l,cr/m,n,cr + /iT,,n+l,cr^m,n,cr) ^ 
m,n m,nja 



where 



I -i6'(m)^ ^* 

^m,n — l,cr ^m,n,cr' 

(6) 



and ^m,n,a = {fL,n,a)- 

When Harm(^, ^) = 0, the Hamiltonian is invariant with re- 
spect to discrete translations along the ?/-axis and the system 
adopts a strip geometry, namely fm,n,cr = fm,a- Since the 
flux $ = p/q is rational, one may set fm,a = fm-\-q,a and the 
2D system reduces to a g x 1 supercell. In order to study the 
effects induced by a harmonic trap, we consider that the har- 
monic potential has its minimum at the centre of the supercell 
and is expressed as 



Harm(^) 



(7) 



inside the supercell, where V2 is the potential strength. 
The mean-field Hamiltonian is given by a sum of single-site 
terms Hmf = X^m ^mp which are represented by 4 x 4 
matrices in the occupation number basis \Nm,i;Nm,i), where 



N„ 



0, 1. The non-zero matrix elements are 



(o,o|w;jp|o,o) = c™, 

(1, 0|W^p|l,0) = (0, 1|HSp|0, 1) = C„ - m + Harm(m), 
(1, 1\H^^\1, 1) = C„ - 2/X + 2Karm(m) + U, 
(0,0|W;Jf|1,0> = -t(*„»+i,T + *m-i,T + 2Cos0(m)«'«,t), 
(0, 0\n^p\0, 1) = -t + + 2 Cos ^(m)*™,x) , 

(o,o|w;Jf|i,o) = {o,i|7i^Sp|i,i), 

(0,0|W;jp|0,l> = {l,0|?iSp|l,l), (8) 



where 



2Cos(9(m)|^, 



m,cr I / 7 



(9) 



and ^(m) = 27r^m. 

The physical quantities of interest for this work are the local 
density pm {Gm\nm\Gm), its variance (n^) - 

{fim)'^ and the average density p = (1/g) where we 



have introduced the notation = J2a '^rn^a- We note that 
Pm since there are two spin states. The ground state \Gm) 
is written as 

\Gm)= O^^'^^'lNm^^NmA). (10) 

where am = {aln'^\ a^'^\ aln'^\ aln'^^ ) is the eigenvector 
corresponding to the lowest eigenvalue of the matrix Eq.([8]). 
The average values = {Gm\f}n,a\Gm) are given by 

*™.i = («rT"r^ + («^^'^T"^;l'°^- (11) 

Eqs. ^ , ([T0| and (TT) form a set of self-consistent equations, 
which are solved in order to evaluate the ground state, the 
average values ^rn,cr and the physical quantities p, pm and 
A^. 

The local density pm is expressed as 

\Gm) 

= |a^)p + |a^)p + 2|aW)|2^ (12) 
the local density's variance is given by 

A™ = {G^\nl,\Gm) - pL 

+ + 2|a2'^f)' (13) 

and the average density p is evaluated according to 

P = - Y{Gm\rim\Gm) 
^ m 

= ^El^^''^l' + l^^n'''^l' + 2|a(^'^^ (14) 



III. NUMERICAL INVESTIGATION 

A. Phase diagram in the /i — t plane for the homogeneous 
system 

In the absence of the harmonic trap V2 = 0, the system is 
invariant under discrete translations along the x direction. As 
a result, the local density Pm = P and its variance = A 
are constant in the supercell. In this case, the representation 
of the average density p in the t — /i plane shows two 
distinctive plateaus (see Fig|2]). The latter correspond to 
the Mott-insulator (resp. band-insulator) phases which are 
characterized by integer fillings, p = I (resp. p = 2), and 
by a vanishing variance A = 0. Contrary to the predictions 
obtained on the basis of quantum Monte-Carlo simulations 
for the ID fermionic optical lattice |8|, the variance A 
evaluated within our mean-field theory plays the role of an 
order parameter for the metal-insulator transition. In the right 
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Figure 2: (Color online) Average density p — p (t, /i) in the t — /x 
plane, for V2 = $ = 0. Mott-insulator and band-insulator phases, 
corresponding to integer fillings, depict two distinctive plateaus at 
p — 1 and p — 2. In the right corner: average density p (blue) and 
average variance A (red) as a function of the hopping amplitude t, 
for /i = 0.7. Inside the Mott lobe p — 1 and A = 0. The hopping 
parameter t and the chemical potential p are expressed in units of the 
on-site interaction energy U . 

corner of Fig [2] we represent the density p and the variance 
A as a function of the hopping parameter t: inside the Mott 
lobe (t < 0.05 for /i = 0.7), p = 1 and the variance indeed 
vanishes A = 0. Outside of this Mott-insulator domain, the 
variance has a finite value which increases as the hopping 
amplitude intensifies and the average density p > 1. 
In order to justify the mean-field approximation Eq. (|4]), we 
have evaluated the fluctuations \Sm,a\ = \{fm,a) - {fm,a)'^\, 
and we have verified that \Sm,a\ <^ {fm,a) for an important 
region surrounding the Mott lobes. It is worth noticing that 
{fm,a) is not an order parameter for the MMIT since it takes 
a non-zero value in a small region within the borders of the 
Mott lobes. 



B. In the presence of a harmonic trap 

When the harmonic trap is added within the supercell m G 
[l^q], the particles agglomerate at the center of the trap and 
the local density is not necessarily constant. In order to study 
this effect, we show the local density pm as a function of the 
chemical potential p in Fig [3] We note that for sufficiently 
large supercells (q > 50), the results obtained in this section 
aren't affected by the specific value given to the parameter q. 
Therefore we choose an arbitrary value and we set q = 153, 
the computations are thus performed in a 153 x 1 supercell, 
and we fix the hopping parameter t = 0.02. In accordance 
with the results obtained with quantum Monte-Carlo simu- 
lations for the ID fermionic optical lattice |8|, we observe 
the emergence of local Mott-insulator phases in the center 
of the system as the chemical potential increases. In Fig [3] 
which corresponds to the case ¥2 = 2, one observes two dis- 




Figure 3: (Color online) Local density pm = Pm (p) as a function of 
the chemical potential p, for t = 0.02 [/, V2 = 2 [/ and for ^ = 153 
lattice sites in the supercell. The chemical potential p is expressed in 
units of the on-site interaction energy U. 

tinctive plateaus of Mott-insulator and band-insulator phases. 
For small values of the chemical potential (p < 0.098 for 
V2 = 2), we observe a parabolic density profile which points 
out that the whole system is in a metallic phase |8, 9|. When 
the chemical potential varies between 0.098 < p < 0.91, 
the system is characterized by a Mott-insulator phase with 
Prn = 1, localized at the center of the trap and surrounded by 
two metallic phases. For 0.91 < p < 1.08, a localized metal- 
lic phase develops in the middle of the Mott-insulator phase. 
When the chemical potential is greater than p > 1.08, a band- 
insulator phase with p^ = 2 appears at the center of the sys- 
tem: the latter is then characterized by three distinctive insu- 
lator phases, which are generally separated and surrounded by 
metallic phases. We point out that such a behaviour also oc- 
curs in the one-dimensional case: Fig[3]reproduces the density 
profiles obtained with quantum Monte Carlo simulations (see 
Fig. 1 and Fig. 2 of Ref. 18J), and with the Bethe ansatz (see 
Fig. 1 of Ref. 13). 



In order to identify the local insulator phases illustrated in 
Fig [3] we represent the local density pm and the variance A^ 
for two fixed values of the chemical potential (p = 0.7 and 
p = 1.7) in Fig|4] In this figure, one observes that a Mott 
plateau indeed develops at the center of the trap for /i = 0.7 
(see Fig|4](a)). This plateau is characterized by prn = 1 and 
by a vanishing variance A^ = 0. We notice that the vari- 
ance A^ is non-zero in the metallic region, where pm is not 
an integer. When p = 1.7, the central plateau corresponds to 
Pm = 2, while two plateaus at p^n = 1 remain near the edges 
(see Fig|4] (6)). The metallic regions, which surround each 
insulator plateau is characterized by a non-vanishing variance 

We also notice that the density vanishes at the edges of the 
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Figure 4: (Color online) Local density pm (blue asterisk) and the 
local variance Am (red dot) when V2 — 2U\ (a) p — 0.7 U and (6) 
/i = 1.7 [/. We set t = 0.02 [/ and ^ = 153 in this figure. 



supercell for /i < 1.8, according to the fact that the parti- 
cles concentrate at the center of the harmonic trap. In fact, 
as the chemical potential increases, the edges undergo differ- 
ent phases: for 1.8 < fi < 2.1 the edges are in a metallic 
phase, for 2.1 < ji < 2.8 the edges are in a Mott-insulator 
phase with pm = 1, for 2.8 < p < 3.1 the edges are again 
in a metallic phase, and for /i > 3.1 the whole system is in a 
global band-insulator phase with p^ = 2. 



As the strength of the harmonic potential V2 increases, the 
global phase diagram p = p{t^p), represented in Fig|2]for 
V2 = 0, varies: the first Mott lobe (p = 1) progressively 
diminishes and eventually disappears for V2 > I, while the 
second insulator region (p = 2) appears at increasing val- 
ues of the chemical potential. This effect is illustrated in 
Figjsj which shows the average density P = \ Xlm Pm 
the average variance A = ^ as a function of the 

chemical potential, for various values of the potential strength 
V2 ^ [O5 2]. Fig|5]is actually a cut through the /i — t phase 
diagram, at a fixed hopping parameter t = 0.02. In this fig- 
ure, one observes that for V2 < I, the system reaches the 
two insulator domain as the chemical potential increases. The 
Mott-insulator and band-insulator phases are characterized by 
a vanishing average variance A. For sufficiently small poten- 
tial strength, the function p = p^p) has two plateaus at p = 1 
and p = 2. The first plateau progressively diminishes while 
V2 is increased (see Fig|5]a, Fig|5]6 and Fig|5]c, correspond- 
ing to V2 = 0, = 0.4 and V2 = 0.8 respectively). When 
V2 > 1, the function p = p{p) only reaches one plateau at 



Figure 5: (Color online) Average density p = ^J2mP^ (^^^ 
ted line) and average variance A = ^J2m^^ (black line) as a 
function of the chemical potential p, for t = 0.02 U and q = 153: 
(a) V2 = 0, (b) V2 = 0.4 U, (c) V2 = 0.8 U, (d) V2 = lU, 
(e) V2 = 1.2U and (/) V2 = 2U. Global Mott-insulator (resp. 
band-insulator) phases correspond to integer fillings p = 1 (resp. 
p = 2) and vanishing variance A = 0. The chemical potential p is 
expressed in units of the on-site interaction energy U. 



p = 2, for /i > 2 (see Fig|5] d, Fig|5] e and Fig|5] /, corre- 
sponding toV2 = l,V2 = 1.2 and ¥2=2 respectively). 
It is worth noticing that a system with non- vanishing average 
variance A may exhibit local insulator phases inside the trap. 
This is the case for the specific situation illustrated in Fig|4] 
(b) which exhibits two plateaus at p^n = 1 and pm = 2, while 
p ^ 1.5 and A ^ 0.02 (for p = 1.7). We finally note the in- 
teresting fact that if the average variance is vanishing, A = 0, 
then the local density pm is constant (pm gives the integer fill- 
ing of the whole lattice) in spite of the inhomogeneous con- 
fining potential. We can conclude that a global Mott-insulator 
phase enhances a homogenization of the system. 



C. Mott metal-insulator phase transition in the presence of a 
"magnetic" field 

The mean-field theory developed in the previous section al- 
lows one to investigate the metal-insulator phase transition in 
a wide range of the parameters. Besides, one of the great ad- 
vantages of our theory resides in the fact that we can investi- 
gate the MMIT when the two-dimensional system rotates. In 
particular, the "high magnetic field" regime ^ ^ 1, can be 
considered. 
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Figure 6: (Color online) Average variance A= ^J^^^Amasa 
function of the hopping parameter t and the chemical potential /i, 
for V2 = 0: (a) $ = and (6) $ = 0.2. The hopping parameter 
t and the chemical potential /x are expressed in units of the on-site 
interaction energy U. 



The MMIT phase diagram, represented in the /i — t plane, is 
affected by the presence of a high "magnetic" field. In particu- 
lar, the Mott lobe corresponding to p = 1 extends over greater 
values of the hopping parameter t. We show this phenomenon 
in Fig|6j which illustrates the phase diagram for ^ = and 
$ = 0.2 in terms of the average variance A (the Mott lobes 
are designated by A = 0). For <l> = 0, the Mott-lobe ex- 
tends to the critical value tc ^ 0.06, whereas it reaches the 
critical value tc ~ 0.08 for <l> = 0.2. We also note that the 
band-insulator domain (p = 2) is modified by the "magnetic" 
field: for /i = 1.3, one observes a critical value tc ^ 0.07 for 
$ = 0, whereas one finds tc ^ 0.1 for <l> = 0.2. Finally we 
point out that the insulator lobes are only modified according 
to the value of the "magnetic" flux ^ = p/q, and that they are 
not affected by the supercell length q. 



A richer representation of this effect is drawn in Fig |7] which 
shows the critical hopping parameter tc as a function of the 
flux ^ = p/q for 3. fixed value of the chemical potential (we 
set /i = 0.3). Here we set q = 227, so that p and q are mu- 
tally primes for all p = 1, 2, ...g 1411 . The phase boundary 
tc = tc{^) is highly irregular and exhibits pronounced slopes 
at ^ = 1/2, $ = 1/3 and <l> = 2/5. This figure is symmetric 
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Figure 7: (Color online) Critical hopping parameter tc as a function 
of the effective "magnetic" flux $ = p/q, for q = 277 and the 
chemical potential /a = 0.3 U. In this figure, we set V2 = 0. The 
hopping parameter t is expressed in units of the on-site interaction 
energy U. 



with respect to $ = 0.5 and is periodic with period = 1. 
The general shape of the phase boundary tc = tc{^) is similar 
for all the values of the chemical potential: in particular, the 
MMIT which takes place around the band-insulator domain 
(p = 2) is also described by such a curve. 
In fact, this important numerical result reveals that the phase 
boundary for the MMIT is directly related to the edge of 
the Hofstadter butterfly represented in Fig[T] The fact that 
this phase transition is intimately connected to the underlying 
single-particle physics is not too surprising: in the context of 
rotating Bose gases, it is already known that the phase bound- 
ary tc = tc{^) for the superfluid-Mott insulator transition also 
displays the Hofstadter butterfly's edge ||28l|34l. We empha- 
size here that our numerical results suggest that this fact holds 
for rotating fermions. This important result enforces the idea 
that rotations, or artificial "magnetic" fields, induces a mod- 
ification to the MMIT in a way which is independent of the 
particles statistics. Moreover the universal character of this 
modification is dictated by the Hofstadter butterfly. As al- 
ready discussed in Ref. |28|, it should be noted that most of 
the atoms are static in the vicinity of the insulator domains. 
The gas of mobile particles and holes is thus highly diluted in 
these regions. This remark suggests that the system behaves 
according to single-particle physics around the insulator re- 
gions and therefore partially explains why these phase bound- 
aries may be related to the fractal single-particle spectrum. 
Finally, we note that when the system reaches the metallic 
phase, A 7^ 0, local insulator phases may still be present 
in the system because of the presence of a "magnetic" field. 
Contrary to the effects induced by a trapping potential and 
described in the previous section, these local Mott-insulator 
phases are generally periodically distributed in space. 
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IV. CONCLUSION 

In this work, we have studied an important feature of 
ultracold Fermions physics. We have shown that two- 
dimensional fermionic optical lattices should exhibit a 
Mott metal-insulator phase transition. In the presence of a 
trapping potential, we have observed the formation of local 
Mott-insulator phases, generally situated at the center of the 
trap. These local Mott-insulator phases are characterized by 
a vanishing variance of the local density. When the system 
is subjected to an external "magnetic" field, or equally when 
the system rotates, we have put forward the fact that the 



Mott lobes are sensible to the specific value of the effective 
"magnetic" flux: for a fixed chemical potential, the critical 
hopping parameter tc = tc{^) indeed follows a fractal 
curve as the "magnetic" flux is varied. As in the case of 
the superfluid-insulator transition in 2D rotating bosonic 
systems, this fractal phase boundary is directly related to the 
single-particle energy spectrum. 
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